/*
 * Open Source Physics software is free software as described near the bottom of this code file.
 *
 * For additional information and documentation on Open Source Physics please see:
 * <http://www.opensourcephysics.org/>
 */

/**
 * particles start in a left box and
 * eventually fill the simulation cell showing the approach to equilibrium. There
 * is a button which allows you to reverse the velocities and see the particles go
 * back to their original positions.
 * @author Jan Tobochnik
 * @author Peter Sibley
 */
package org.opensourcephysics.stp.neto;

import org.opensourcephysics.display.Drawable;
import org.opensourcephysics.display.DrawingPanel;

import java.awt.*;
import java.util.HashMap;
import java.util.Map;

public class Diamagnetism implements Drawable {
    public double x[], y[];
    public Vector v[],v0[];
    public Vector a[];
    public Map<Integer, ParticulaColidindo> colidindo;  
    public Color color[];    
    public int numberOfParticles = 1;
    public double cellLength; //tamanho da caixa
    public double rho = 0.2; // fixed density
    public double x0;
    public double temperature = 5.0;
    public int steps = 0;
    public double dt = 0.01;
    public int n0 = 0, n1 = 0;
    double radius = 0.3;
    public double t = 0;
    public double b = 5; //magnetic field
    public double m = 1; //mass of electron
    public double q = 1; //charge of electron
    public double totalMCentro = 0; //magnetizacao total no centro
    public double totalMBorda = 0; //magnetizacao total na borda
    

    public Vector vecB; //magnetic field

    /**
     * Constructor Diamagnetism         momento angular grift cap 5 ou 7 bor eo b mi =  fator 2 mi*0b/r0 fator 2
     */
    public Diamagnetism() {
        setArrays();
    }

    /**
     * Metodo principal que faz os calculos a cada passo
     * */
    public void step() { // Velocity Verlet algorithm
        double dt2half = 0.5 * dt * dt;
        n0 = 0;
        n1 = 0;
        for (int i = 0; i < numberOfParticles; i++) {
            x[i] += v[i].x * dt + a[i].x * dt2half;
            y[i] += v[i].y * dt + a[i].y * dt2half;
            x[i] = image(x[i]);
            y[i] = image(y[i]);
            v[i].x += a[i].x * dt;
            v[i].y += a[i].y * dt;
            //normaliza
            v[i] = v[i].normalize().produtoEscalar(v0[i].modulo());
        }
        wallColision();
        accel();
        calcMagnetizacao();

        steps++;
        t += dt;      
    }

    private void calcMagnetizacao() {
        //To change body of created methods use File | Settings | File Templates.
        totalMCentro = 0;
        totalMBorda = 0;
        //Faz uma aproximacao desse circulo imaginario usando a igualdade do periodo
        double raioBorda = (4*(cellLength-radius))/(Math.PI*2);
        for (int i = 0; i < numberOfParticles; i++) {
            ParticulaColidindo p = colidindo.get(i);

            if (p==null){
                //esta girando no centro
                //calcula o raio
                double r = m / q * v[i].modulo() / b;
                //magnetizacao
                totalMCentro += -0.5*q*v[i].modulo()*r;
            } else {
                //esta batendo nas bordas
                totalMBorda += 0.5*q*(p.distanciaTotal/t)*raioBorda; //na direcao de z 
            }
        }
        //System.out.println("totalMCentro=" + totalMCentro + " totalMBorda=" + totalMBorda + " total=" + (totalMCentro+totalMBorda));
    }
    
    private void accel() {
        for (int i = 0; i < numberOfParticles; i++) {
            Vector f = v[i].produtoVetorial(vecB).produtoEscalar(q);
            //f=m.a --> a = f/m
            a[i] = f.produtoEscalar(1/m);
        }
    }


    /**
     * Isso deve ser periodic boundary conditions, utilizado para determinar a colisao
     */
    private double pbc(double s) {
        if (s > 0.5 * cellLength) {
            s -= cellLength;
        } else if (s < -0.5 * cellLength) {
            s += cellLength;
        }
        return s;
    }


    /**
     * Esse metodo que cria a parede infinita
     */
    private double image(double s) {
        if (s > cellLength) {
            s -= cellLength;
        } else if (s < 0) {
            s += cellLength;
        }
        return s;
    }

    private void wallColision() {
        for (int i = 0; i < numberOfParticles; i++) {
            //verifica se o x atingiu alguma das paredes
            if (x[i] + radius > cellLength || x[i] - radius < 0) {
                v[i].x = -1 * v[i].x;
                calculaDistanciaColicao(i);
                color[i] = Color.RED;
            }
            //verifica se o y atingiu alguma das paredes
            if (y[i] + radius > cellLength || y[i] - radius < 0) {
                v[i].y = -1 * v[i].y;
                calculaDistanciaColicao(i);
                color[i] = Color.RED;
            }
        }
    }

    private void calculaDistanciaColicao(int i) {
        if (colidindo.get(i)==null){
            colidindo.put(i, new ParticulaColidindo(x[i],y[i]));
        } else {
            colidindo.get(i).calculaAdicionaDistancia(x[i],y[i]);
        }
    }

    public void reverse() {
        for (int i = 0; i < numberOfParticles; i++) {
            v[i].x = -v[i].x;
            v[i].y = -v[i].y;            
        }
    }

    public void draw(DrawingPanel myWorld, Graphics g) {
        if (x == null) {
            return;
        }
        int pxRadius = Math.abs(myWorld.xToPix(radius) - myWorld.xToPix(0));
        int pyRadius = Math.abs(myWorld.yToPix(radius) - myWorld.yToPix(0));
        if (pxRadius < 1) {
            pxRadius = 1;
        }
        if (pyRadius < 1) {
            pyRadius = 1;
        }
        for (int i = 0; i < numberOfParticles; i++) {
            int xpix = myWorld.xToPix(x[i]) - pxRadius;
            int ypix = myWorld.yToPix(y[i]) - pyRadius;
            g.setColor(color[i]);
//      Color color = Color.BLACK;
//      if(x[i]>cellLength/2) color=Color.BLUE;
//      g.setColor(color);
            g.fillOval(xpix, ypix, 2 * pxRadius, 2 * pyRadius);
        }
        g.setColor(Color.black);
        int xpix = myWorld.xToPix(0);
        int ypix = myWorld.yToPix(0);
        int lx = myWorld.xToPix(cellLength) - myWorld.xToPix(0);
        int ly = myWorld.yToPix(cellLength) - myWorld.yToPix(0);
        g.drawRect(xpix, ypix, lx, ly);
        g.drawLine(xpix, ypix, xpix, ypix + ly);
        g.drawLine(xpix + lx, ypix, xpix + lx, ypix + ly);
        g.drawLine(xpix, ypix, xpix + lx, ypix);
        g.drawLine(xpix, ypix + ly, xpix + lx, ypix + ly);
    }


    // ---------------- Inicializacao ----------------------------------------------------------------------------------
    
    public void setArrays() {
        int N = numberOfParticles;
        x = new double[N];
        y = new double[N];
        v = new Vector[N];
        v0 = new Vector[N];
        a = new Vector[N]; 
        color = new Color[N];
        colidindo = new HashMap<Integer, ParticulaColidindo>();
    }

    public void initialize() {
        vecB = new Vector(0,0,b);
        cellLength = Math.sqrt(numberOfParticles / rho);
        System.out.println("cellLength=" + cellLength);
        reset();
        setArrays();
        setPositions();
        setVelocities();
        setColor();
        accel();
    }

    public void setPositions() {
        boolean overlap;
        double ds = cellLength; //preenche a caixa inteira
        for (int i = 0; i < numberOfParticles; ++i) {
            do {
                overlap = false;
                do{
                    x[i] = ds * Math.random();
                } while (x[i] < radius || x[i] > cellLength - radius); //garante que a bola nao fica presa na borda

                do{
                    y[i] = cellLength * Math.random();
                } while (y[i] < radius || y[i] > cellLength - radius); //garante que a bola nao fica presa na borda

                int j = 0;
                while ((j < i) && !overlap) { //isso evita que as bolas ja comecem colidindo
                    double dx = x[i] - x[j];
                    double dy = y[i] - y[j];
                    dy = pbc(dy);
                    if (dx * dx + dy * dy < 1.0) {
                        overlap = true;
                    }
                    j++;
                }
            } while (overlap);
        }
    }

    private void setColor() {
        for (int i = 0; i < numberOfParticles; ++i) {
            //color[i] = new Color((float) Math.random() ,(float) Math.random(),(float) Math.random());
            color[i] = new Color((float)0.0, (float)0.5, (float)0.0); //green
        }
    }

    public void setVelocities() {
        double tpi = 2.0 * Math.PI;
        for (int i = 0; i < numberOfParticles; ++i) {
            v0[i] = new Vector();
            double r1 = Math.random();
            double r2 = Math.random() * tpi;
            v0[i].x = Math.sqrt(-2 * temperature * Math.log(r1)) * Math.cos(r2); // sets random initial velocities
            v0[i].y = Math.sqrt(-2 * temperature * Math.log(r1)) * Math.sin(r2);
        }
        // zero Momentum
        double vxsum = 0.0;
        double vysum = 0.0;
        for (int i = 0; i < numberOfParticles; ++i) {
            vxsum += v0[i].x;
            vysum += v0[i].y;
        }
        double vxcm = vxsum / numberOfParticles;
        double vycm = vysum / numberOfParticles;
        for (int i = 0; i < numberOfParticles; ++i) {
            v0[i].x -= vxcm;
            v0[i].y -= vycm;

            v[i] = new Vector();
            v[i].x = v0[i].x; 
            v[i].y = v0[i].y;
        }
    }   

    public void reset() {
        steps = 0;
        t = 0;
        n0 = 0;
        n1 = 0;
        x0 = cellLength / 2.0;
    }

    public void zeroAverages() {
        steps = 0;
        n0 = 0;
        n1 = 0;
    }
}


class Vector {
    public double x = 0;
    public double y = 0;
    public double z = 0;

    public Vector() {

    }

    public Vector(double x, double y, double z) {
        this.x = x;
        this.y = y;
        this.z = z;
    }

    public double modulo() {
        return Math.sqrt(x * x + y * y + z * z);
    }

    public Vector normalize() {
        double r = modulo();
        return new Vector(x / r, y / r, z / r);
    }

    public Vector produtoVetorial(Vector v) {
        return new Vector(y * v.z - z * v.y, z * v.x - x * v.z, x * v.y - y * v.x);
    }

    public Vector produtoEscalar(double v) {
        return new Vector(x * v, y * v, z * v);
    }
}

class ParticulaColidindo{
    public double lastYColisao = 0;
    public double lastXColisao = 0;
    public double distanciaTotal = 0;

    ParticulaColidindo(double lastXColisao, double lastYColisao) {
        this.lastYColisao = lastYColisao;
        this.lastXColisao = lastXColisao;
    }

    public void calculaAdicionaDistancia(double currentXColisao, double currentYColisao){
        distanciaTotal += Math.sqrt( Math.pow(lastXColisao-currentXColisao,2) + Math.pow(lastYColisao-currentYColisao,2));
        lastXColisao = currentXColisao;
        lastYColisao = currentYColisao;
    }
}

/*
 * Open Source Physics software is free software; you can redistribute
 * it and/or modify it under the terms of the GNU General Public License (GPL) as
 * published by the Free Software Foundation; either version 2 of the License,
 * or(at your option) any later version.

 * Code that uses any portion of the code in the org.opensourcephysics package
 * or any subpackage (subdirectory) of this package must must also be be released
 * under the GNU GPL license.
 *
 * This software is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston MA 02111-1307 USA
 * or view the license online at http://www.gnu.org/copyleft/gpl.html
 *
 * Copyright (c) 2007  The Open Source Physics project
 *                     http://www.opensourcephysics.org
 */
